Juan M. Fonseca-Solís · Agosto 2020 · 5 min read
Los estetoscopios digitales son herramientas usadas por los médicos para realizar remotamente un seguimiento de sus pacientes por problemas de salud relacionados, por ejemplo, con el corazón y los pulmones. En muchos de los casos, estos estetoscopios digitales utilizan el códec OPUS, el mismo códec empleado en sistemas de teleconferencia como Zoom, MS Teams y Skype. El éxito del OPUS y otros códecs similares radica en el hecho de que son capaces de transmitir voz y música de alta calidad empleando bajas tasas de bits; para ello eliminan ciertas frecuencias que el oído humano es incapaz de escuchar. Estas frecuencias eliminadas podrían ser necesarias para los algoritmos de reconocimiento de patrones. En este ipython notebook realizamos pruebas al códec OPUS para determinar si la pérdida de la calidad es significativa. Los resultados obtenidos son preliminares, pero muestran que aplicando un filtrado pasabajas y usando una tasa de bits de 8 kbps la distorsión es mínima.
Figura 0. Tomada de https://stethoscopeshop.eu/16833-thickbox_default/thinklabs-one-digital-stethoscope.jpg.
El códec OPUS fue propuesto en el 2013 durante la conferencia número 135 de la Audio Engineering Society, como una alternativa más eficiente respecto de otros algoritmos como el MP3, Advance Audio Coding (AAC) y Vorbis. La figura 1 muestra que el códec OPUS ofrece la misma calidad de audio que estos algoritmos pero a una menor tasa de bits, lo que permite ofrecer calidades mejores con el mismo ancho de banda, y es capaz de cubrir casi todo el rango de transmisión. OPUS está compuesto por dos algoritmos llamados SILK (creado por Skype Limited para comprimir voz) y Constrained Energy Lapped Transform (CELT) (para comprimir música) [3, 4, 13]; así como un modo híbrido.
Figura 1. Tomada de https://opus-codec.org/static/comparison/quality.svg.
Antes de continuar con el diseño de los casos de prueba para evaluar al OPUS, es necesario explicar algunos concepto de psicoacústica. En primer lugar, es importante mencionar que códecs como el OPUS son algoritmos de compresión "con pérdida", esto quiere decir que reducen el tamaño de las grabaciones eliminando información auditiva que de por sí es inaudible para el oído humano. La forma de identificar estas frecuencias inaudibles es estudiando la cóclea, y entendiendo que este órgano sufre dos efectos llamados enmascaramiento frecuencial y enmascaramiento temporal. Para efectos de este ipython solo nos interesa el primero. Como muestra la figura 2, la cóclea es un órgano en forma cónica encargado de percibir las frecuencias y transmitirlas como impulsos eléctricos al cerebro. Se encuentra "arrollado" dentro del oído interno y recibe estímulos del oído externo e intermedio mediante la ventana oval (esta se conecta con los huesecillos conectados al tímpano).
Figura 2. Tomada de https://en.wikipedia.org/wiki/Auditory_system#/media/File:Anatomy_of_the_Human_Ear.svg.
El enmascaramiento frecuencial es un concepto relacionado con unas escalas llamadas Bark y ERB, que definen las denominadas bandas críticas. Las escalas Bark y ERB son modelos matemáticos utilizados para representar la resolución frecuencial de la cóclea. Se basan en el hecho de que el oído humano realiza un mejor trabajo reconociendo tonos graves que tonos agudos (distribución logarítmica). La escala Equivalent Rectangular Bandwidth (ERB) es una versión mejorada de la escala Bark y su ecuación es $\text{ERB}(f)=21.4 \log_{10}{(0.00437f+1)}$, mostrada a continuación [9].
%pylab inline
import numpy as np
def erb(F):
return 21.4 * np.log10(0.00437*F+1)
F = np.linspace(0,20000,1000) # 0 a 20k Hz es el rango de audición humana
E = erb(F)
pylab.plot(F/1e3,E)
pylab.ylabel('Frequency (kHz)')
pylab.xlabel('ERB (kHz)')
pylab.show()
Por su parte, las bandas críticas son segmentos de la cóclea, de diferente longitud y localización, en los cuales el oído humano no es capaz de diferenciar dos tonos lo suficientemente cercanos. Dos ejemplos de bandas críticas se ubica en el rango [555,677] kHz con frecuencia central de 511 Hz y [8928,10353] Hz con frecuencia central 9589 Hz [10].
E = np.linspace(0,erb(20000),32) # 32 puntos equidistantes
F = (np.power(10,E/21.4)-1)/0.00437 # invertimos el logaritmo de la escala ERB
print(F.astype(int))
pylab.stem(F,np.ones(len(F)))
pylab.xlabel('Frecuencia central de cada banda crítica (Hz)')
pylab.show()
La figura 3 ilustra que el enmascaramiento frecuencial se produce cuando en una misma banda crítica la frecuencia de mayor energía enmascara a las de menor energía usando su envolvente espectral.
Figura 3. Tomada de https://output.com/blog/9-sound-design-tips-to-hack-your-listeners-ears.
El objetivo propuesto de verificar que el códec OPUS sea apto para un algoritmo de reconocimiento de patrones no se logra verificando que la compresión del audio se haga sin pérdida (pues la teoría psicoacústica de la sección anterior nos dice que el oido ya realiza ese truncamiento de información) sino más bien verificando que la pérdida de información no sea tan brusca que deje al espectro irreconocible. Considerando que el rango de frecuencias deseado es de $[60, 1000]$ Hz para los pulmones y $[20, 500]$ Hz para el corazón ($[20, 1000]$ Hz en conjunto), luego de haber comprimido y descomprimido el audio, podemos definir los siguientes casos de prueba [14, 15]:
Para probar el caso de prueba #1 se puede utilizar un barrido de frecuencias como sigue [6]: $$ F(t) = \Big(\frac{F_1-F_0}{T}\Big)t + F_0, $$
donde el barrido puede ser construido usando un senoidal de la forma $x[t] = A \sin(2\pi Ft)$. El rango de frecuencias a usar es de $[20, 20000]$ Hz (el rango de audición humana).
%pylab inline
from scipy.io import wavfile
from IPython.display import Audio
import numpy as np
def plotSpecgram(x,Fs,fMax=None):
# graficar el espectrograma
fig, ax = pylab.subplots(nrows=1)
ax.specgram(x, NFFT=1024, Fs=Fs, noverlap=900)
pylab.xlabel('Tiempo (s)')
pylab.ylabel('Frecuencia (Hz)')
if fMax != None:
pylab.ylim([0,fMax])
pylab.show()
def plotFFT(x,Fs,xlim=None):
pylab.figure()
N2=int(len(x)/2)
F = np.linspace(0,Fs/2,N2)/1e3
X = np.sqrt(np.abs(np.fft.fft(x)[0:N2]))
pylab.plot(F, X)
X[0] = 0 # remover el valor 0 que contiene el promedio de la energía de la magnitud espectral
pylab.xlabel('Frecuencia (kHz)')
pylab.ylabel('$\sqrt{|S(f)|}$')
if xlim != None:
pylab.xlim(xlim)
pylab.show()
rango = [20.0, 20000.0] # el rango promedio de audicion humano en Hz
Fs = 44.1 * 1e3 # la tasa de muestreo de los equipos comerciales
T = 1.0 # segundos (t1-t0)
N = int(T*Fs)
n = np.arange(0,N)
F0 = (rango[1]-rango[0])*n/N + rango[0]
x = np.sin(np.pi*F0/Fs*n) # f=F0/Fs: frecuencia discreta
plotSpecgram(x,Fs)
plotFFT(x,Fs)
wavfile.write('/tmp/barrido.wav',int(Fs),x.astype(np.float32))
Audio(x, rate=Fs)
El espectrograma muestra que el barrido va de los 20 a los 20000 Hz como se deseaba. Asimismo, la magnitud espectral muestra una recta horizontal casi en todo el rango, indicando que todas las frecuencias están presentes.
Procedemos a codificar y decodificar el barrido de frecuencias usando OPUS para una tasa de bits de 8 Kbps (la misma empleada en una videollamada entre dos personas) [12].
def readPlayVisualizeFile(inputFile,fMax=None):
fs, x = wavfile.read(inputFile)
y = np.array(x)/max(x)
if None==fMax:
plotSpecgram(x,fs)
else:
plotSpecgram(x,fs,fMax)
return fs, x
#!sudo apt-get install ffmpeg
#!sudo apt-get install opus-tools
!ffmpeg -loglevel error -y -i /tmp/barrido.wav -qscale 0 /tmp/wavRaw.wav
!opusenc --quiet --bitrate 8 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav
Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)
Audio(x, rate=Fs)
El espectrograma muestra que partir de los casi 5 kHz el barrido de frecuencia se pliega, en un efecto llamado aliasing, lo que produce que veamos frecuencias que en realidad no existen, esto sí que podría afectar la calidad del clasificador. Para resolverlo se podría emplear un filtro pasabajas con frecuencia de corte 5 kHz antes de realizar la compresión. Por su parte, la magnitud espectral muestra que a partir de los 5 kHz las frecuencias altas son atenuadas. Probemos ahora con 128 kbps (fullband stereo) a ver si la calidad mejora. Sin embargo, aunque se introdujeron artefactos que antes no estaban, el rango $[20,1000]$ Hz sigue estando casi intacto, lo cual indica que para un algoritmo de reconocimiento de patrones lo suficientemente selectivo, la corrupción en la señal no debería afecta (hay que considerar si existen armónicas para asegurar esto completamente).
Procedemos a codificar y decodificar el barrido de frecuencias usando OPUS para una tasa de bits de 128 Kbps (full-band) para determinar si con una tasa de bits más alta la calidad mejora [12].
!ffmpeg -loglevel error -y -i /tmp/barrido.wav -qscale 0 /tmp/wavRaw.wav # opus-tools
!opusenc --quiet --bitrate 128 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav
Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)
Audio(x, rate=Fs)
En efecto, la calidad ha mejorado. El espectrograma sigue mostrando que hay cierto aliasing pero con bastante menor energía, lo cual favorable; y la magnitud espectral muestra una recta horizontal similar a la que vimos con la señal cruda. Por lo que podemos decir que con full-band ofrece una menor distorsión.
Nota: En este punto valdría la pena preguntarse cuál de los tres métodos de OPUS realizó la compresión y si los resultados variarían en caso de que se fuerce algún método en particular (SILK, CELT o el híbrido).
Para el caso de prueba 2, podemos construir una señal de dos tonos como sigue:
y = 0.1*np.sin(np.pi*595/(Fs/2)*n) + np.sin(np.pi*605/(Fs/2)*n)
plotSpecgram(y,Fs)
plotFFT(y,Fs,xlim=[500/1e3,700/1e3])
wavfile.write('/tmp/enmascaramiento.wav',int(Fs),y.astype(np.float32))
Audio(y, rate=Fs)
Ahora procedemos a realizar la compresión y decompresión a una tasa de 8 kbps:
!ffmpeg -loglevel error -y -i /tmp/enmascaramiento.wav -qscale 0 /tmp/wavRaw.wav # opus-tools
!opusenc --quiet --bitrate 8 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav
Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs,xlim=[500/1e3,700/1e3])
Audio(x, rate=Fs)
El espectrograma muestra que hay información nueva a partir de los 5 kHz que no contenía la grabación original. Esto no parece ser producto de un aliasing sino de información agregada por el OPUS, y que podría eliminarse usando un filtrado pasa bajas post-procesamiento. Se observa que el tono de 605 Hz sigue siendo dominante y que el de 595 Hz se ha enmascarado levemente, pero que aún sigue siendo distinguible.
Procedemos ahora a codificar y decodificar una grabación real tomada de la Respiratory Sound Database usando OPUS con la misma tasa de bits de 8 Kbps [16]. La grabación corresponde a un paciente femenina de 75 años ausculada en la traquea con un micrófono de grabación esterea de modelo AKG C417L.
Nota: a partir de este momento se recomienda usar audifonos para apreciar mejor la calidad el audio.
Fs, x = readPlayVisualizeFile('./wav/107_3p2_Tc_mc_AKGC417L_2.wav')
plotFFT(x,Fs)
Audio(x, rate=Fs)
Observamos que hay energía alta en el rango $[0,1.5]$ y otra menor en $[7.5,9.5]$ kHz. La primera parece corresponder a los sonidos naturales del cuerpo, la segundo para más bien una frecuencia parásita que podría haberse colado al momento de realizar la grabación (el cuerpo humano no produce tonos constantes y puros como los que se observan después de los 7500 Hz).
# OPUS
!ffmpeg -loglevel error -y -i ./wav/107_3p2_Tc_mc_AKGC417L_2.wav -qscale 0 /tmp/wavRaw.wav # same quality
!opusenc --quiet --bitrate 8 /tmp/wavRaw.wav /tmp/opusEnc.opus
!opusdec --quiet /tmp/opusEnc.opus /tmp/opusDec.wav
Fs, x = readPlayVisualizeFile('/tmp/opusDec.wav')
plotFFT(x,Fs)
Audio(x, rate=Fs)
El espectrograma muestra que la energía un poco antes de los 5 kHz fue atenuada, pero ya descubrimos que esto no importa tanto pues se sale del rango objetivo en los $[20,1000]$ Hz. Lo más importante parece ser que las frecuencias están inalterada en el rango deseado. La magnitud espectral muestra un espectro limpio.
Los resultados arrojados por los casos de prueba nos permiten decir que la calidad del OPUS es aceptable para la tasa de bits usada en las llamadas 1:1 (de 8kps) cuando los sonidos a analizar corresponden al corazón y los pulmones. Antes de emplear OPUS con una tasa de bits de 8kps, es necesario aplicar un filtro pasabajas con frecuencia de corte 5 kHz para evitar el efecto del aliasing y evitar también confundir al algoritmo de reconocimiento de patrones con información que no existe. El seguir estas dos consideraciones parece ser una condición suficiente para emplear el OPUS sin tener que llegar a transmitir datos sin comprimir (lo cual sería difícil técnicamente). Otra alternativa sería usar el OPUS a una tasa de 128 kbps donde se encontró corrupción casi nula. Como trabajo futuro, se podría investigar un método para cuantificar la distorsión apreciada en los casos de prueba propuestos y contemplar más casos de prueba (por ejemplo, para el enmascaramiento temporal y evaluar una colección de sonidos).

Este obra está bajo una licencia de Creative Commons Reconocimiento-NoComercial-SinObraDerivada 4.0 Internacional. El sitio juanfonsecasolis.github.io es un blog dedicado a la investigación independiente en temas relacionados al procesamiento digital de señales. Para reutilizar este artículo y citar las fuente por favor utilice el siguiente Bibtex:
@online{Fonseca2020,
author = {Juan M. Fonseca-Solís},
title = { Pruebas por pares o pairwise testing},
year = 2020,
url = {https://juanfonsecasolis.github.io/blog/JFonseca.pairwisetesting.html},
urldate = {}
}